// <copyright file="MergedLines.cs" company="Microsoft Corporation">
// The MIT License (MIT)
// 
// Copyright (c) 2014, Microsoft Corporation
// 
// Permission is hereby granted, free of charge, to any person obtaining a copy
//  of this software and associated documentation files (the "Software"), to deal
//  in the Software without restriction, including without limitation the rights
//  to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
//  copies of the Software, and to permit persons to whom the Software is
//  furnished to do so, subject to the following conditions:
// 
// The above copyright notice and this permission notice shall be included in
//  all copies or substantial portions of the Software.
// 
// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
//  IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
//  FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
//  AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
//  LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
//  OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN
//  THE SOFTWARE.
// </copyright>
// --------------------------------------------------------------------------------------------------------------------

namespace Microsoft.Robotics.Vision.ImageProcessing2D
{
    using System;
    using System.Collections.Generic;

    using Microsoft.Robotics.Numerics;
    using Microsoft.Robotics.Tracing;

    /// <summary>
    /// A merged line is a merge of other merged lines or image lines segments
    /// Merging is rule driven. Lines will merge if they satisfy
    /// 1) Have same slope
    /// 2) Endpoints are in close proximity
    /// 3) The length of the merged line is close to the sum of length of the segments
    /// </summary>
    public class MergedLines
    {
        /// <summary>
        /// Unique Identifier for line
        /// </summary>
        public readonly int Idx;

        /// <summary>
        /// Angle threshold in degrees
        /// </summary>
        private static double angleMergeThreshold = 18;

        /// <summary>
        /// Distance fraction threshold between endpoints for merge consideration
        /// </summary>
        private static double lengthMergeThreshold = 80;

        /// <summary>
        /// Distance fraction threshold between endpoints for merge consideration
        /// </summary>
        private static double endPointMergeThreshold = 0.05;

        /// <summary>
        /// Tracing context for merged lines
        /// </summary>
        private static TraceContext mergeTraceContext = TraceContexts.FromString(TraceContexts.Vision, "ImageLineMerge");

        /// <summary>
        /// Counter for assigning line index
        /// </summary>
        private static int indexCounter = 0;

        /// <summary>
        /// Initializes a new instance of the <see cref="MergedLines" /> class
        /// </summary>
        /// <param name="line">First element in the merge</param>
        public MergedLines(ImageLine line)
        {
            this.Start = new Point2D<double>(line.StartX, line.StartY);
            this.End = new Point2D<double>(line.EndX, line.EndY);
            this.AverageEdgeStrength = line.AverageEdgeStrength;
            this.VarEdgeStrength = line.VarEdgeStrength;
            this.MeanEdgelDeviation = line.MeanEdgelDeviation;
            this.VarEdgelDeviation = line.VarEdgelDeviation;
            this.Components = new List<MergedLines>();
            this.TotalLength = MergedLines.LineEndDistance(this.Start.X, this.Start.Y, this.End.X, this.End.Y);
            this.Idx = System.Threading.Interlocked.Increment(ref MergedLines.indexCounter);
        }

        /// <summary>
        /// Initializes a new instance of the <see cref="MergedLines" /> class
        /// </summary>
        /// <param name="mergedLine">Initial merged line</param>
        public MergedLines(MergedLines mergedLine)
        {
            this.Start = new Point2D<double>(mergedLine.Start.X, mergedLine.Start.Y);
            this.End = new Point2D<double>(mergedLine.End.X, mergedLine.End.Y);
            this.AverageEdgeStrength = mergedLine.AverageEdgeStrength;
            this.VarEdgeStrength = mergedLine.VarEdgeStrength;
            this.MeanEdgelDeviation = mergedLine.MeanEdgelDeviation;
            this.VarEdgelDeviation = mergedLine.VarEdgelDeviation;
            this.Components = new List<MergedLines>();
            this.TotalLength = mergedLine.TotalLength;
            this.IsMerged = mergedLine.IsMerged;
            this.Idx = System.Threading.Interlocked.Increment(ref MergedLines.indexCounter);
        }

        /// <summary>
        /// Gets or sets the angle threshold for considering line merges in degree
        /// </summary>
        public static double AngleMergeThreshold
        {
            get
            {
                return MergedLines.angleMergeThreshold;
            }

            set
            {
                MergedLines.angleMergeThreshold = Math.Max(value, 1);
            }
        }

        /// <summary>
        /// Gets or sets the length fraction threshold for considering merges
        /// </summary>
        public static double LengthMergeThreshold
        {
            get
            {
                return MergedLines.lengthMergeThreshold;
            }

            set
            {
                MergedLines.lengthMergeThreshold = Math.Max(value, 1);
            }
        }

        /// <summary>
        /// Gets or sets the distance fraction threshold between endpoints for merge consideration
        /// </summary>
        public static double EndPointMergeThreshold
        {
            get
            {
                return MergedLines.endPointMergeThreshold;
            }

            set
            {
                MergedLines.endPointMergeThreshold = Math.Min(Math.Max(value, 0), 1);
            }
        }

        /// <summary>
        /// Gets Component lines
        /// </summary>
        public List<MergedLines> Components { get; private set; }

        /// <summary>
        /// Gets angle in degrees line makes with y axis in range [0-180]
        /// </summary>
        public double Angle
        {
            get
            {
                return MergedLines.GetLineAngle(this.Start, this.End);
            }
        }

        /// <summary>
        /// Gets binned distance of line start from origin
        /// </summary>
        public int BinnedStartDistance
        {
            get
            {
                return (int)(MergedLines.LineEndDistance(this.Start.X, this.Start.Y, 0, 0) / MergedLines.LengthMergeThreshold);
            }
        }

        /// <summary>
        /// Gets binned distance of line End from Origin
        /// </summary>
        public int BinnedEndDistance
        {
            get
            {
                return (int)(MergedLines.LineEndDistance(this.End.X, this.End.Y, 0, 0) / MergedLines.LengthMergeThreshold);
            }
        }

        /// <summary>
        /// Gets the binned angle using bin truncation
        /// </summary>
        public int BinnedAngle
        {
            get
            {
                return MergedLines.FloatAngleToBin(this.Angle, 0.0);
            }
        }

        /// <summary>
        /// Gets the binned angle with bin rounding up
        /// </summary>
        public int BinnedAngleRoundUp
        {
            get
            {
                return MergedLines.FloatAngleToBin(this.Angle, 1.0);
            }
        }

        /// <summary>
        /// Gets distance of line normal to origin
        /// </summary>
        public double DistanceFromOrigin
        {
            get
            {
                double delX = this.End.X - this.Start.X;
                double delY = this.End.Y - this.Start.Y;
                double n = Math.Sqrt(delX * delX + delY * delY);
                return (delX * this.Start.Y - delY * this.Start.X) / n;
            }
        }

        /// <summary>
        /// Gets Total merged length
        /// </summary>
        public double TotalLength { get; private set; }

        /// <summary>
        /// Gets Start point of merged line
        /// </summary>
        public Point2D<double> Start { get; private set; }

        /// <summary>
        /// Gets End point of merged line
        /// </summary>
        public Point2D<double> End { get; private set; }

        /// <summary>
        /// Gets the average edgel strength
        /// </summary>
        public double AverageEdgeStrength { get; private set; }

        /// <summary>
        /// Gets the Variance of edgel strength
        /// </summary>
        public double VarEdgeStrength { get; private set; }

        /// <summary>
        /// Gets the mean deviation between the edge data and the line fit
        /// </summary>
        public double MeanEdgelDeviation { get; private set; }

        /// <summary>
        /// Gets the variance of deviation between data and the
        /// fitted line
        /// </summary>
        public double VarEdgelDeviation { get; private set; }

        /// <summary>
        /// Gets or sets a value indicating whether the line has merged
        /// </summary>
        public bool IsMerged { get; set; }

        /// <summary>
        /// Euclidean distance between two points
        /// </summary>
        /// <param name="p1x">Point 1 X</param>
        /// <param name="p1y">Point 1 Y</param>
        /// <param name="p2x">Point 2 X</param>
        /// <param name="p2y">Point 2 Y</param>
        /// <returns>Euclidean Distance</returns>
        public static double LineEndDistance(double p1x, double p1y, double p2x, double p2y)
        {
            double x = p1x - p2x;
            double y = p1y - p2y;
            return Math.Sqrt(x * x + y * y);
        }

        /// <summary>
        /// Get line angle from line end points
        /// </summary>
        /// <param name="start">Line start</param>
        /// <param name="end">Line end</param>
        /// <returns>Angle in degrees</returns>
        public static double GetLineAngle(Point2D<double> start, Point2D<double> end)
        {
            double delX = end.X - start.X;
            double delY = end.Y - start.Y;
            double angle = Math.Atan2(delY, delX);
            angle = angle < 0.0 ? angle + Math.PI : angle;
            return angle / MathConstants.Degrees2Radians;
        }

        /// <summary>
        /// Bin a floating point angle to binned angle using current binning parameters
        /// </summary>
        /// <param name="angle">Floating point angle to convert</param>
        /// <param name="roundFactor">Rounding factor</param>
        /// <returns>Binned angle</returns>
        public static int FloatAngleToBin(double angle, double roundFactor)
        {
            int binAngle = (int)(angle / MergedLines.AngleMergeThreshold + roundFactor);
            int maxAngles = (int)((180 / MergedLines.AngleMergeThreshold));
            binAngle = binAngle % maxAngles;
            return binAngle;
        }

        /// <summary>
        /// Runs line merging on a set of 2D input line segments obtained from an image. Returns a set containing (without duplication) 
        /// greedy permutations of line merging that satisfy merging criteria. Merging criteria consist 
        /// of threshold on distance and orientation matches between candidate segments. 
        /// As an example if there are 3 input line segments that can be merged the returned set contains
        /// the original set plus (line1 + line2), (line1 + line2 + line3)
        /// NOTE: this method is not thread safe.
        /// </summary>
        /// <param name="lines">Lines to merge</param>
        /// <param name="angleMergeThreshold">Merge angle threshold (degree)</param>
        /// <param name="normMergeThreshold">Normalized merging distance threshold, expressed as a fraction of merging line candidates</param>
        /// <param name="endPointMergeThreshold">Threshold distance between line segment end points when considering merges</param>
        /// <param name="minLength">Minimum input line segment length to consider</param>
        /// <param name="maxBinAngles">Number of binned angles when doing merging. Essentially sets an angle threshold for merging</param>
        /// <returns>List of merged lines</returns>
        public static List<MergedLines> DoLineMerge(List<ImageLine> lines, double angleMergeThreshold, double normMergeThreshold, double endPointMergeThreshold, double minLength, int maxBinAngles)
        {
            MergedLines.AngleMergeThreshold = angleMergeThreshold;
            MergedLines.LengthMergeThreshold = normMergeThreshold;
            MergedLines.EndPointMergeThreshold = endPointMergeThreshold;

            Dictionary<Bin2D, List<MergedLines>> hist = new Dictionary<Bin2D, List<MergedLines>>();
            List<MergedLines> mergedLines = new List<MergedLines>();

            lines.Sort((a, b) => { return -a.Length.CompareTo(b.Length); });

            // For first round of merging create a histogram
            // by line segment angle and end distance from origin
            foreach (ImageLine l in lines)
            {
                if (l.Length < minLength)
                {
                    continue;
                }

                MergedLines lm = new MergedLines(l);
                mergedLines.Add(lm);

                MergedLines.AddToHist(hist, lm.BinnedAngle, lm.BinnedStartDistance, lm);
                if (lm.BinnedEndDistance != lm.BinnedStartDistance)
                {
                    MergedLines.AddToHist(hist, lm.BinnedAngle, lm.BinnedEndDistance, lm);
                }

                if (lm.BinnedAngleRoundUp != lm.BinnedAngle)
                {
                    MergedLines.AddToHist(hist, lm.BinnedAngleRoundUp, lm.BinnedStartDistance, lm);
                    if (lm.BinnedEndDistance != lm.BinnedStartDistance)
                    {
                        MergedLines.AddToHist(hist, lm.BinnedAngleRoundUp, lm.BinnedEndDistance, lm);
                    }
                }
            }

            int mergeCount = 0;

            // Keep merging until no merge occur
            // The first round uses binning by angle and end point distance. Subsequent 
            // rounds only bin by angle
            for (int iter = 0; mergeCount > 0 || iter < 2; ++iter)
            {
                MergedLines.DumpHist(hist);
                mergeCount = MergedLines.DoMergeRound(hist, endPointMergeThreshold, mergedLines);
                if (mergeCount > 0 || iter < 1)
                {
                    hist.Clear();
                    HashSet<int> added = new HashSet<int>();

                    // Rebin all lines
                    foreach (MergedLines lm in mergedLines)
                    {
                        if (false == lm.IsMerged && false == added.Contains(lm.Idx))
                        {
                            MergedLines.AddToHist(hist, lm.BinnedAngle, 0, lm);
                            if (lm.BinnedAngle != lm.BinnedAngleRoundUp)
                            {
                                MergedLines.AddToHist(hist, lm.BinnedAngleRoundUp, 0, lm);
                            }

                            added.Add(lm.Idx);
                        }
                    }
                }
            }

            return mergedLines;
        }

        /// <summary>
        /// Checks if this line intersects another line
        /// </summary>
        /// <param name="other">The other line</param>
        /// <param name="maxDistance">Max tolerance on intersection distance</param>
        /// <returns>True if lines intersect</returns>
        public bool DoLinesIntersect(MergedLines other, double maxDistance)
        {
            bool doIntersect = (MergedLines.LineEndDistance(this.Start.X, this.Start.Y, other.Start.X, other.Start.Y) < maxDistance) ||
                                         (MergedLines.LineEndDistance(this.End.X, this.End.Y, other.Start.X, other.Start.Y) < maxDistance) ||
                                         (MergedLines.LineEndDistance(this.Start.X, this.Start.Y, other.End.X, other.End.Y) < maxDistance) ||
                                         (MergedLines.LineEndDistance(this.End.X, this.End.Y, other.End.X, other.End.Y) < maxDistance);

            return doIntersect;
        }

        /// <summary>
        /// Executes a single round of line merging. Merging is a greedy operation
        /// </summary>
        /// <param name="hist">Histogram of lines</param>
        /// <param name="mergeDistanceRelThreshold">Distance merging threshold expressed as fraction of line length</param>
        /// <param name="mergedLines">Set of merged lines. Newly merged lines are added to the list</param>
        /// <returns>Number of lines merged</returns>
        private static int DoMergeRound(Dictionary<Bin2D, List<MergedLines>> hist, double mergeDistanceRelThreshold, List<MergedLines> mergedLines)
        {
            int mergeCount = 0;

            ////TraceOut.Info(MergedLines.mergeTraceContext, "\nMerging\nAngBin\tDistBin\tAngle\tIdx\tStartX\tStartY\tEndX\tEndY\tLength\tStrength\tAngle\tIdx\tStartX\tStartY\tEndX\tEndY\tLength\tStrength");

            HashSet<int> mergedSet = new HashSet<int>();

            foreach (KeyValuePair<Bin2D, List<MergedLines>> kvp in hist)
            {
                List<MergedLines> linesInBin = hist[kvp.Key];

                for (int i = 0; i < linesInBin.Count; ++i)
                {
                    MergedLines ml = linesInBin[i];
                    if (true == mergedSet.Contains(ml.Idx))
                    {
                        continue;
                    }

                    for (int j = 0; j < linesInBin.Count; ++j)
                    {
                        MergedLines line = linesInBin[j];
                        if (mergedSet.Contains(line.Idx))
                        {
                            continue;
                        }

                        MergedLines merged = ml.DoesMerge(line, mergeDistanceRelThreshold);
                        if (null != merged)
                        {
                            line.IsMerged = true;
                            ml.IsMerged = true;
                            mergedLines.Add(merged);
                            ////TraceOut.Info(MergedLines.mergeTraceContext, "{0}\t{1}\t{2:F1}\t{3}\t{4:F1}\t{5:F1}\t{6:F1}\t{7:F1}\t{8:F1}\t{9:F1}\t{10:F1}\t{11}\t{12:F1}\t{13:F1}\t{14:F1}\t{15:F1}\t{16:F1}\t{17:F1}", kvp.Key.Bin1, kvp.Key.Bin2, ml.Angle, ml.Idx, ml.Start.X, ml.Start.Y, ml.End.X, ml.End.Y, ml.TotalLength, ml.AverageEdgeStrength, line.Angle, line.Idx, line.Start.X, line.Start.Y, line.End.X, line.End.Y, line.TotalLength, line.AverageEdgeStrength);

                            if (i > j)
                            {
                                --i;
                            }

                            linesInBin.RemoveAt(j);
                            mergedSet.Add(line.Idx);
                            mergedSet.Add(ml.Idx);
                            ++mergeCount;
                            break;
                        }
                    }
                }
            }

            return mergeCount;
        }

        /// <summary>
        /// Add a candidate line to the binned histogram
        /// </summary>
        /// <param name="hist">Parent histogram</param>
        /// <param name="binnedAngle">Binned angle key</param>
        /// <param name="distBin">Binned distance from origin key</param>
        /// <param name="lineEx">Line to add</param>
        private static void AddToHist(Dictionary<Bin2D, List<MergedLines>> hist, int binnedAngle, int distBin, MergedLines lineEx)
        {
            Bin2D ls = new Bin2D(binnedAngle, distBin);

            if (true == hist.ContainsKey(ls))
            {
                List<MergedLines> linesInBin = hist[ls];
                linesInBin.Add(lineEx);
            }
            else
            {
                hist.Add(ls, new List<MergedLines>() { lineEx });
            }
        }

        /// <summary>
        /// Debug dump of merged line histogram
        /// </summary>
        /// <param name="hist">Histogram to dump</param>
        private static void DumpHist(Dictionary<Bin2D, List<MergedLines>> hist)
        {
            // Debug info
            ////TraceOut.Info(MergedLines.mergeTraceContext, "\nIdx\tBinAngle\tAngle\tDistOri\tLen\tStartX\tStartY\tEndx\tEndY\tstrength");
            foreach (KeyValuePair<Bin2D, List<MergedLines>> kvp in hist)
            {
                List<MergedLines> linesInBin = hist[kvp.Key];

                for (int i = 0; i < linesInBin.Count; ++i)
                {
                    MergedLines l = linesInBin[i];
                    ////TraceOut.Info(MergedLines.mergeTraceContext, "{0}\t{1}\t{2:F1}\t{3}\t{4:F1}\t{5:F1}\t{6:F1}\t{7:F1}\t{8:F1}\t{9:F1}", l.Idx, kvp.Key.Bin1, l.Angle, kvp.Key.Bin2, l.TotalLength, l.Start.X, l.Start.Y, l.End.X, l.End.Y, l.AverageEdgeStrength);
                }
            }
        }

        /// <summary>
        /// Checks if a candidate merges with current line. Merging criteria include 
        /// - Distance between start / end points between candidate and current line
        /// - Size of gap if merging occurs
        /// - Change in line slope if merging occurs
        /// </summary>
        /// <param name="line">Candidate line</param>
        /// <param name="mergeDistanceRelThreshold">Merging threshold expressed for endpoints of line expressed as a fraction of current line length</param>
        /// <returns>The new merged line or null if no merge occurred</returns>
        private MergedLines DoesMerge(MergedLines line, double mergeDistanceRelThreshold)
        {
            bool didMerge = false;

            Point2D<double> start = new Point2D<double>(this.Start);
            Point2D<double> end = new Point2D<double>(this.End);

            double dist = MergedLines.LineEndDistance(this.Start.X, this.Start.Y, line.Start.X, line.Start.Y);
            double newLength = this.TotalLength + line.TotalLength;
            double maxOldLen = Math.Max(this.TotalLength, line.TotalLength);
            double minDistance = newLength * mergeDistanceRelThreshold;

            if (dist < minDistance)
            {
                double endPointLen = MergedLines.LineEndDistance(this.End.X, this.End.Y, line.End.X, line.End.Y);
                if (Math.Abs(endPointLen - newLength) < minDistance && endPointLen > maxOldLen)
                {
                    didMerge = true;
                    start = new Point2D<double>(line.End.X, line.End.Y);
                    minDistance = dist;
                }
            }

            dist = MergedLines.LineEndDistance(this.Start.X, this.Start.Y, line.End.X, line.End.Y);
            if (dist < minDistance)
            {
                double endPointLen = MergedLines.LineEndDistance(this.End.X, this.End.Y, line.Start.X, line.Start.Y);
                if (Math.Abs(endPointLen - newLength) < minDistance && endPointLen > maxOldLen)
                {
                    didMerge = true;
                    start = new Point2D<double>(line.Start.X, line.Start.Y);
                    end = new Point2D<double>(this.End);
                    minDistance = dist;
                }
            }

            dist = MergedLines.LineEndDistance(this.End.X, this.End.Y, line.Start.X, line.Start.Y);
            if (dist < minDistance)
            {
                double endPointLen = MergedLines.LineEndDistance(this.Start.X, this.Start.Y, line.End.X, line.End.Y);
                if (Math.Abs(endPointLen - newLength) < minDistance && endPointLen > maxOldLen)
                {
                    didMerge = true;
                    end = new Point2D<double>(line.End.X, line.End.Y);
                    start = new Point2D<double>(this.Start);
                    minDistance = dist;
                }
            }

            dist = MergedLines.LineEndDistance(this.End.X, this.End.Y, line.End.X, line.End.Y);
            if (dist < minDistance)
            {
                double endPointLen = MergedLines.LineEndDistance(this.Start.X, this.Start.Y, line.Start.X, line.Start.Y);
                if (Math.Abs(endPointLen - newLength) < minDistance && endPointLen > maxOldLen)
                {
                    didMerge = true;
                    end = new Point2D<double>(line.Start.X, line.Start.Y);
                    start = new Point2D<double>(this.Start);
                    minDistance = dist;
                }
            }

            if (true == didMerge)
            {
                // Check the angle change
                double newAngle = MergedLines.GetLineAngle(start, end);
                double diffAngle = newAngle - this.Angle;

                // Currently not using line segment gradient differences but would like to keep the option
                // of including them
                ////double strengthStat = this.DifferenceMeanTest(this.AverageEdgeStrength, this.VarEdgeStrength, this.TotalLength, line.AverageEdgeStrength, line.VarEdgeStrength, line.TotalLength);
                ////double deviationStat = this.DifferenceMeanTest(this.MeanEdgelDeviation, this.VarEdgelDeviation, this.TotalLength, line.MeanEdgelDeviation, line.VarEdgelDeviation, line.TotalLength);

                if ((Math.Abs(diffAngle) < MergedLines.AngleMergeThreshold ||
                    Math.Abs(diffAngle - 180) < MergedLines.AngleMergeThreshold))
                {
                    ////Math.Abs(strengthStat) < 2.0 &&
                    ////Math.Abs(deviationStat) < 2.0)
                    MergedLines merged = new MergedLines(this);

                    merged.UpdateWithMerge(start, end, line);
                    return merged;
                }
                else
                {
                    didMerge = false;
                }
            }

            return null;
        }

        /// <summary>
        /// Updates the parameters of a newly merged line
        /// </summary>
        /// <param name="start">New start point</param>
        /// <param name="end">New end point</param>
        /// <param name="other">Other line that it merged with</param>
        private void UpdateWithMerge(Point2D<double> start, Point2D<double> end, MergedLines other)
        {
            this.Start = start;
            this.End = end;
            this.Components.Add(other);
            this.TotalLength = MergedLines.LineEndDistance(this.Start.X, this.Start.Y, this.End.X, this.End.Y);
            double varNew;
            double meanNew;
            this.UpdateVarianceMean(this.AverageEdgeStrength, this.VarEdgeStrength, this.TotalLength, other.AverageEdgeStrength, other.VarEdgeStrength, other.TotalLength, out varNew, out meanNew);
            this.AverageEdgeStrength = meanNew;
            this.VarEdgelDeviation = varNew;
            this.UpdateVarianceMean(this.MeanEdgelDeviation, this.VarEdgelDeviation, this.TotalLength, other.MeanEdgelDeviation, other.VarEdgelDeviation, other.TotalLength, out varNew, out meanNew);
            this.MeanEdgelDeviation = meanNew;
            this.VarEdgelDeviation = varNew;
        }

        /// <summary>
        /// Update mean and variance as sum of two components
        /// </summary>
        /// <param name="v1">First component variance</param>
        /// <param name="m1">First component mean</param>
        /// <param name="n1">First component count</param>
        /// <param name="v2">Second component variance</param>
        /// <param name="m2">Second component mean</param>
        /// <param name="n2">Second component count</param>
        /// <param name="v">Variance of sum</param>
        /// <param name="m">Mean of sum</param>
        private void UpdateVarianceMean(double v1, double m1, double n1, double v2, double m2, double n2, out double v, out double m)
        {
            v = 0;
            m = 0;
            double n = n1 + n2;

            if (n <= 0)
            {
                return;
            }

            m = ((m1 * n1) + (m2 * n2)) / n;

            v = ((v1 + (m1 * m1)) * n1 + (v2 + (m2 * m2) * n2)) / n - m * m;
        }

        /// <summary>
        /// Simple z tests of difference between means
        /// </summary>
        /// <param name="m1">Mean of sample 1</param>
        /// <param name="v1">Variance of sample 1</param>
        /// <param name="n1">Sample 1 count</param>
        /// <param name="m2">Mean of sample 2</param>
        /// <param name="v2">Variance of sample 2</param>
        /// <param name="n2">Sample 2 count</param>
        /// <returns>Test statistic</returns>
        private double DifferenceMeanTest(double m1, double v1, double n1, double m2, double v2, double n2)
        {
            double sampVar = Math.Sqrt((v1 / n1) + (v2 / n2));
            double z = (m1 - m2) / sampVar;

            return z;
        }

        /// <summary>
        /// Index into 2D sparse array
        /// </summary>
        internal struct Bin2D
        {
            /// <summary>
            /// Line angle
            /// </summary>
            public readonly int Bin1;

            /// <summary>
            /// Line norm
            /// </summary>
            public readonly int Bin2;

            /// <summary>
            /// Initializes a new instance of the <see cref="Bin2D" /> struct
            /// </summary>
            /// <param name="bin1">First bin value</param>
            /// <param name="bin2">Second bin value</param>
            public Bin2D(int bin1, int bin2)
            {
                this.Bin1 = bin1;
                this.Bin2 = bin2;
            }
        }
    }
}
